{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# [StatsLib](https://github.com/kthohr/stats)\n",
    "\n",
    "StatsLib is a templated C++ library of statistical distribution functions, featuring unique compile-time computing capabilities and seamless integration with several popular linear algebra libraries.\n",
    "\n",
    "Features:\n",
    "* A header-only library of probability density functions, cumulative distribution functions, quantile functions, and random sampling methods.\n",
    "* Functions are written in a specialized C++11 `constexpr` format, enabling the library to operate as both a compile-time and run-time computation engine.\n",
    "* Designed with a simple **R**-like syntax.\n",
    "* Optional vector-matrix functionality with wrappers to support:\n",
    "    * [Armadillo](http://arma.sourceforge.net/)\n",
    "    * [Blaze](https://bitbucket.org/blaze-lib/blaze)\n",
    "    * [Eigen](http://eigen.tuxfamily.org/index.php)\n",
    "* Matrix-based operations are parallelizable with OpenMP.\n",
    "* Released under a permissive, non-GPL license.\n",
    "\n",
    "Author: [Keith O'Hara](https://www.kthohr.com)\n",
    "\n",
    "License: Apache Version 2\n",
    "\n",
    "## Distributions\n",
    "\n",
    "Functions to compute the cdf, pdf, quantile, as well as random sampling methods, are available for the following distributions:\n",
    "\n",
    "* Bernoulli\n",
    "* Beta\n",
    "* Binomial\n",
    "* Cauchy\n",
    "* Chi-squared\n",
    "* Exponential\n",
    "* F\n",
    "* Gamma\n",
    "* Inverse-Gamma\n",
    "* Laplace\n",
    "* Logistic\n",
    "* Log-Normal\n",
    "* Normal (Gaussian)\n",
    "* Poisson\n",
    "* Student's t\n",
    "* Uniform\n",
    "* Weibull\n",
    "\n",
    "In addition, pdf and random sampling functions are available for several multivariate distributions:\n",
    "\n",
    "* inverse-Wishart\n",
    "* Multivariate Normal\n",
    "* Wishart\n",
    "\n",
    "## This Notebook\n",
    "\n",
    "To run a code cell, use `shift + enter`.\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 1,
   "metadata": {},
   "outputs": [],
   "source": [
    "// include libraries\n",
    "#include <iostream>              // for printing\n",
    "#include \"../include/stats.hpp\""
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Syntax and Examples\n",
    "\n",
    "Functions are called using an **R**-like syntax. Some general rules:\n",
    "\n",
    "* density functions: `stats::d*`. For example, the Normal (Gaussian) density is called using\n",
    "``` cpp\n",
    "stats::dnorm(<value>,<mean parameter>,<standard deviation>);\n",
    "```\n",
    "* cumulative distribution functions: `stats::p*`. For example, the Gamma CDF is called using\n",
    "``` cpp\n",
    "stats::pgamma(<value>,<shape parameter>,<scale parameter>);\n",
    "```\n",
    "* quantile functions: `stats::q*`. For example, the Beta quantile is called using\n",
    "``` cpp\n",
    "stats::qbeta(<value>,<a parameter>,<b parameter>);\n",
    "```\n",
    "* random sampling: `stats::r*`. For example, to generate a single draw from the Logistic distribution:\n",
    "``` cpp\n",
    "stats::rlogis(<location parameter>,<scale parameter>,<seed value or random number engine>);\n",
    "```\n",
    "\n",
    "Some examples:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "metadata": {},
   "outputs": [],
   "source": [
    "// evaluate the normal PDF at x = 1, mu = 0, sigma = 1\n",
    "double dval_1 = stats::dnorm(1.0,0.0,1.0);\n",
    " \n",
    "// evaluate the normal PDF at x = 1, mu = 0, sigma = 1, and return the log value\n",
    "double dval_2 = stats::dnorm(1.0,0.0,1.0,true);\n",
    " \n",
    "// evaluate the normal CDF at x = 1, mu = 0, sigma = 1\n",
    "double pval = stats::pnorm(1.0,0.0,1.0);\n",
    " \n",
    "// evaluate the Laplacian quantile at p = 0.1, mu = 0, sigma = 1\n",
    "double qval = stats::qlaplace(0.1,0.0,1.0);\n",
    "\n",
    "// draw from a t-distribution dof = 30\n",
    "double rval = stats::rt(30);"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "0.241971\n",
      "-1.41894\n",
      "0.841345\n",
      "-1.60944\n",
      "0.20657\n"
     ]
    }
   ],
   "source": [
    "std::cout << dval_1 << std::endl;\n",
    "std::cout << dval_2 << std::endl;\n",
    "std::cout << pval << std::endl;\n",
    "std::cout << qval << std::endl;\n",
    "std::cout << rval << std::endl;"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Seeding\n",
    "\n",
    "Random number seeding is available in two formats: seed values and random number engines.\n",
    "\n",
    "* Seed values are passed as unsigned integers. For example, to generate a draw from a normal distribution N(1,2) with seed value 1776:\n",
    "``` cpp\n",
    "stats::rnorm(1,2,1776);\n",
    "```\n",
    "* Random engines in StatsLib use the 64-bit Mersenne-Twister generator (`std::mt19937_64`) and are passed by reference. Example:\n",
    "``` cpp\n",
    "std::mt19937_64 engine(1776);\n",
    "stats::rnorm(1,2,engine);\n",
    "```"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "metadata": {
    "scrolled": true
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "random draws: 3.17135, 3.17135\n"
     ]
    }
   ],
   "source": [
    "stats::ullint_t seed_val = 1776UL;\n",
    "\n",
    "double ran_val_1 = stats::rnorm(1,2,seed_val);\n",
    "\n",
    "std::mt19937_64 engine(seed_val);\n",
    "double ran_val_2 = stats::rnorm(1,2,engine);\n",
    "\n",
    "std::cout << \"random draws: \" << ran_val_1 << \", \" << ran_val_2 << std::endl;"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Compile-time Computation Capabilities\n",
    "\n",
    "StatsLib is designed to operate equally well as a compile-time computation engine. Compile-time computation allows the compiler to replace function calls (e.g., `dnorm(0,0,1)`) with static values in the source code. That is, functions are evaluated during the compilation process, rather than at run-time. This capability is made possible due to the templated `constexpr` design of the library and can be verified by inspecting the assembly code generated by the compiler. \n",
    "\n",
    "The compile-time features are enabled using the `constexpr` specifier. The example below computes the pdf, cdf, and quantile function of the Laplace distribution:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 5,
   "metadata": {},
   "outputs": [],
   "source": [
    "constexpr double cdens  = stats::dlaplace(1.0,1.0,2.0); // answer = 0.25\n",
    "constexpr double cprob  = stats::plaplace(1.0,1.0,2.0); // answer = 0.5\n",
    "constexpr double cquant = stats::qlaplace(0.1,1.0,2.0); // answer = -2.218875..."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 6,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "result: 0.25, 0.5, -2.21888\n"
     ]
    }
   ],
   "source": [
    "std::cout << \"result: \"<< cdens << \", \" << cprob << \", \" << cquant << std::endl;"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "C++11",
   "language": "C++11",
   "name": "xeus-cling-cpp11"
  },
  "language_info": {
   "codemirror_mode": "text/x-c++src",
   "file_extension": ".cpp",
   "mimetype": "text/x-c++src",
   "name": "c++",
   "version": "-std=c++11"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 2
}
